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Abstract 

Combining recent moment and sparse semidefinite programming (SDP) relaxation 
techniques, we propose an approach to find smooth approximations for solutions of prob- 
lems involving nonlinear differential equations. Given a system of nonlinear differential 
equations, we apply a technique based on finite differences and sparse SDP relaxations 
for polynomial optimization problems (POP) to obtain a discrete approximation of its 
solution. In a second step we apply maximum entropy estimation (using moments of 
a Borel measure associated with the discrete solution) to obtain a smooth closed-form 
approximation. The approach is illustrated on a variety of linear and nonlinear ordinary 
differential equations (ODE), partial differential equations (PDE) and optimal control 
problems (OOP), and preliminary numerical results are reported. 
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1 Introduction 



Problems involving nonlinear differential equations arise in a variety of models for real world 
problems. Even finding approximate solutions for nonlinear differential equations remains a 
challenge. In the previous work Mevissen et al.(2008), Mevissen et al.(2009)| we have estab- 



lished a technique based on sparse semidefinite programming (SDP) relaxations to construct 
discrete approximations for solutions of systems of nonlinear differential equations. In this 
paper we present a novel approach to obtain smooth approximations for solutions of dif- 
ferential equations and of problems involving differential equations such as optimal control 
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problems. Namely, an approximate solution is obtained by applying the maximum entropy 
estimation of Borwein et al.(1991), Lasserre(2007)| with a finite number of moments of a 



Borel measure associated with a discrete approximation of the solution of the differential 
equation. For linear differential equations, an SDP relaxation based method was proposed in 
Bertsimas et al. (2006)| to generate contracting sequences of lower and upper bounds for the 



moments. 



Our contribution in this paper is primarily concerned with nonlinear differential equations 
and nonlinear optimal control problems. In a first step, we take advantage of the discrete 
approximations provided by the SDP relaxation method in Mevissen et al.(2008)| to compute 



a finite set of moments for an appropriately defined measure with discrete support. Next, 
an approximation for the solution of the differential equation is obtained in closed form by 
maximum entropy estimation, using the moments of the discrete measure. To the best of our 
knowledge, it seems to be the first attempt to apply maximum entropy estimation to obtain 
smooth approximations for solutions of linear and nonlinear differential equations. Finally, if 
maximum entropy estimation does guarantee some weak convergence of the estimate to the 
true solution as the number of moments increases, it does not guarantee pointwise convergence 
on the entire domain of the differential equation. However, as our preliminary results show 
for different linear and nonlinear differential equations and optimal control problems, accurate 
pointwise approximations can be achieved on certain regions of the domain. 

The structure of the paper is as follows. In Section [2] we briefly recall the sparse SDP relax- 
ations method for solving nonlinear differential equations numerically, which is the basis for 
our technique. In Section [3] we introduce the method of maximum entropy estimation. Our 
technique to compute smooth approximations for the solutions of ordinary differential equa- 
tions (ODE), partial differential equations (PDE) and optimal control problems (OCP), which 
combines the methods from Section [2] and [3l is presented in Section HI Finally, preliminary 
numerical results for this method are reported in Section |5l 



2 Sparse SDP relaxations for solving differential equa- 
tions 



In this section we recall the approach to compute discrete approximations for solutions to 
systems of differential equations with polynomial data presented in Mevissen et al.(2008)| . 



2.1 Transforming a differential equation into a POP 

In this paper we consider bidimensional differential equation problems of the following type 

L{u{x, y)) + G{u{x, y)) = f{x, y) \/{x, y) e il, 

H{u{x,y)) = g{x,y) y{x,y) e dn, (1) 
Ibd < u{x,y) < ubd W{x,y) G Q, 

with Q := [a:min, Xmax] X [l/miin 2/max] C. M^, f : Q ^ M™, g : dQ — 7- smooth functions, 
L : P — 7- "D and H : V ^ V differential operators, G : I? — t- I? operator, L, H and G are 
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polynomial in the function m G P and its derivatives for some function space T> C Li{Vl), with 
Li(f2) the Banach space of integrable functions on f2, and Ibd < ubd bounds for the function 
values of u. Even if lower and upper bounds are not given by the problem directly, we need 
to impose them, since they are crucial for numerical stability of our method. If little is known 
about the solutions of ([1]), loose lower and upper bounds are to be imposed. 

In a first step of the sparse semidefinite programming relaxation (SDPR) method for solving 
differential equations, we transform an ODE or a PDE of form ([T]) into a polynomial opti- 
mization problem (POP). We then discretize the rectangular domain VLhy N (one-dimensional 
case) or N^Ny (two-dimensional case) grid points, approximate the derivates by standard fi- 
nite differences, for instance ^ %^^2^^^ ~ M,+i,j-2tt^,j+ti,-i,j ^ ^ff]^QYe Uij := u{xi,yj), and denote 
the discretized differential operators at some interior grid point {xi,yj) and boundary grid 
point {xkiVi) as Ljj and Hj^i, respectively. Finally, we take the discretization of the problem 
([1]) as a system of constraints, choose some objective function F which needs to be polynomial 
in M, i.e. F{u) G R[u], and obtain the following optimization problem 

min F{u) 

s.t. Lij{u) + G{uij) = f{xi, yj), 

\hdij < Ui,j < ubdij V(i,j) G {1..N,} X {l..Ny}, (2) 
Hk,i{u) = gixk,yi) V(A:, /) G {1,N,} x {L.iVj 

U{1..N,} X {l,Ny} 

which, since all functions are polynomial in the variable u = (mi,i, . . . lUN^^Ny)-, is a POP. In 
the case of an OCP the objective function F is given by the discretization of the optimal value 
function. In the case of a differential equation with several solutions we can pick a specific 
solution by choosing an appropriate objective function. A particular choice for the objective 
function may be motivated by the underlying problem, such as energy functionals in physics. 



2.2 Sparse SDP relaxations 



Problem (121) is of the form: 



POP min F{u) 

s.t. ^,(m)>0 Vj G {!,..., fc}, 

Hu) = Q V2 G {!,...,/}, 

Ibds <Us< ubds Vs G {1, . . . , n}. 



(3) 



where F{u) = XloeN" ^au", gj(u) = XloeN" ^nd hi{u) = J2aen" ^iaU". As mentioned 

above, the dimension of this POP is n = N^Ny (or fi = in the one-dimensional case). 
But due to the structure of the finite difference discretization only a small number of the h 
components of u occurs in each constraint or in each monomial of the objective function of 
(12]). Thus, the POPs derived from nonlinear differential equations are sparse. 



A systematic way of characterizing the structured sparsity of a POP has been introduced in 
[Waki et al.(2006)"] , and this structured sparsity is exploited by the sparse SDP relaxations 
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constructed in 



Waki et al.(2"006)] . For the POP © one obtains the sparse SDP relaxations: 



SDP min Yj\a\ 

s.t. M,^_,^^{gjyJtij)) > V j G {1, . . . , A;} 

M^_^Xh^y,Iii^))=^ €{!,...,[} (4) 
M^{y,lt)>0 \fte{l,...,d} 
Ibds < ?/s < ubds Vs G {1, . . . ,fi}, 

where wj := [^^^|^] , Wi := f^"!^' ] , {It}f=i is the set of subsets of {1, ... , h} derived from the 
correlative sparsity pattern matrix of the POP Waki et al.(2006)], M^(-, ■) and Mw^w.{-, ■) are 



^ ^ ^ ^ ^ _ _ . . ^ " /J ' LU \ 7 / UJ UJj\ J / 

the partial moments and localizing matrices of Lasserre (2006)1 , and 



is the order of the SDP relaxation. As w increases, solving SDP^ generates a nondecreasing 
sequence of lower bounds for min(POP), namely: 

min (SDP^,„,J < min (SDP^^,^+i) < . . . < min (POP) . 

Moreover, if ([3D has a unique minimizer x*, the vector (y* . . . , -) obtained from an optimal 
solution of (jlj) is an approximation of x*, and under certain compactness conditions, see 
Lasserre (2006)1 , as w — )■ CX3, 



min (SDP^) ^ min (POP) and yl x\ (5) 

In implementations of the sparse SDP relaxation such as SparsePOP, small linear perturbation 
terms are added to the objective function F, in order to ensure that ([2]) has a unique optimal 
solution. Therefore, an optimal solution of the SDP (jlj) is an approximation of the optimal 
solution of the POP (EI). 



2.3 Discrete approximation 

From dS]), asymptotic convergence of minimum and minimizer of SDP^ to minimum and min- 
imizer of POP ([3]) are guaranteed under uniqueness of the optimal solution and compactness 
of the feasible set of the POP. However, recall that the dimension n of POP ([3j) is A^^ Ny, and 
so, if w is much larger than lOmin, the resulting SDP^„ in (jl]) is untractable for general purpose 
SDP solvers. Therefore, we have to restrict ourselves to w G {'Wmin,'Wmin + 1} for medium 
scale A^, and we cannot expect that SDP^^.^ or SDP^^^^.^^+i provide accurate approximations 
to a solution of a nonlinear differential equation. In order to improve accuracy, we may apply 
additional locally convergent optimization techniques. For instance sequential quadratic pro- 
gramming (SQP) can be applied to ([2]), starting from the solution of SDP^ as initial guess. 
Since the solution of the SDP relaxation is an approximation for the global optimizer, it is 
a systematic choice for an initial point for local methods which does not require any a priori 
information. Combining the sparse SDP relaxation with SQP (or another local method) is 
summarized in the scheme: 



Method 1 The SDPR method 



1. Choose a discretization (A^, A^) for the differential equation problem. 

2. Choose a relaxation order w and an objective function F. Apply SDP^ to ^ and obtain 
its solution u and the lower bound min (SDP^) for min (POP). 
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3. Apply sequential quadratic programming (SQP) to ^ with u as initial guess, and obtain 
u as discrete approximation to a solution of the differential equation problem. 



3 Maximum entropy estimation 



In this section we briefly introduce the maximum entropy estimation of Borwein et al.(1991) 



Lasserre(200'7)j |Lasserre(2009)| , our second tool to find smooth approximations for solutions 



of nonlinear differential equations. The maximum entropy estimation is concerned with the 
following problem: Let u G Li{VL) be nonnegative and partially known by the finite vector m 
of moments up to order M of the associated Borel measure rf/x := udxdy on VL. From the 
only knowledge of m, find an estimate um G Li[VL) such that all moments of order up to M 
of the measure d^M '■= UMdxdy match those of d^x and analyze the asymptotic behavior of 
Um as M — 7- oo. An elegant method consists in finding the estimate um that maximizes the 
Boltzmann- Shannon entropy functional 



j ulnudxdy, u G -Li(fi). 
Jn 

In this case, the optimal estimate Up^ is given by 



u*M{x,y) =ulj{v*,x,y) = exp ^ v-jx'y^, (6) 



0<i,j<M 
i+j<M 



where v* G M'™' is a (global) optimizer of the convex finite-dimensional optimization problem: 



max(m, t>)— / um{v, x,y) dx dy. (7) 



The optimization problem ([7]) can be solved by first or second order methods like Newton's 
method or SQP. For these methods, gradient and Hessian of the objective function in ([7]) 
must be computed. In the case the domain Q is simple, we may compute them by quadrature 
and cubature formulas. For more difficult domains one may use the procedure described in 



Lasserre(2007) , Bertsimas et al.(2008)j . Concerning the behavior of u\.j as M — )■ oo and its 



relationship with u, one has the following weak convergence result from [Borwein et al.(1991)| 

Proposition 1 Let u\,j be obtained from an optimal solution of ([^. Then, as M ^ oo, 

ijj{x,y) u\,j{x,y) dxdy — > / ip{x,y) u{x,y) dxdy, 
n Jn 

for every bounded measurable function i(j : Q ^ which is continuous almost everywhere. 
However, the pointwise convergence u\j{x,y) — ?■ u{x,y) does not hold in general. 
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4 Smooth approximations to solutions of differential equa- 
tions 



In this section we show how to combine the SDPR method and maximum entropy estimation to 
obtain smooth approximations for solutions of hnear and nonhnear differential equations. Our 
discussion focuses on the 2-dimensional case, but the 1-dimensional case is covered analogously. 
Let u E Vhe a solution of ([T]) and assume without loss of generality that Ibd > so that u is 
nonnegative on Q. For Ibd < define u := u — Ibd and apply the outlined procedure to the 
new function u. Associated with m, let dfi{x,y) := udxdy be the finite Borel measure on Q 
with moment vector m = {rriij) of all moments up to order M: 



rriij := / x^y^d^{x,y)= / x^y^u{x,y)dx dy, 
Jn Jn 

for € with i + j < M. For linear PDEs, a hierarchy of tightening lower and upper 

bounds for the components niij of the moment vector can be obtained by solving a sequence 
of SDPs as proposed in Bertsimas et al. (2006)| . However, this approach cannot be applied in 
the case of nonlinear PDEs. Our strategy is to use a discrete approximation {uk,i)i<k<Nj:,i<i<Ny 
of a solution of the PDF to approximate the moments rrii j by: 



J] J]42/z «A:,;AxA?/, (8) 



A;=l 1=1 



for e W with i + j < M, where Ax := ^""^ , Ay := 2/2|^^aia. If the discretization 

{Nx,Ny) is sufficiently fine and {uk^i)k,i is a close approximation of u, then we expect m^^ to 
be a good approximation of rrii j for all G N^. Thus, we can apply maximum entropy 
estimation to the vector to obtain a smooth approximation of u. This idea is formalized in 
the following algorithm for obtaining smooth approximations of solutions of linear or nonlinear 
PDF problems. 



Method 2 The smooth SDP approximation method 
Given a PDF problem of form ([7]j. 



1. Choose a discretization {Nj.,Ny), relaxation order w and objective F, and apply the 
SDPR method to obtain a discrete approximation {ukj)k,i to a solution of the PDF prob- 
lem. If not given in the formulation of the PDF problem, impose lower and upper bounds, 
Ibd and ubd for u. 

2. Choose a moment bound M G N and use {uk^i)k,i to calculate as in 

3. Apply maximum entropy estimation to and obtain vector v* E RI™^!, optimal solu- 
tion of 

4- Obtain the approximation u\.^ with u\[[x,y) = j exp(t;*ja;*y-') for a solution u of the 
PDF problem (QP on Vt. 
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As M — )■ oo, -u^j — > u weakly (see Proposition [T]) but not pointwise, i.e., one cannot guarantee 
u\,f{x,y) — )■ u{x,y) on fi. Nevertheless, as reported in [Lasserre(2007)| the maximum entropy 
estimation may provide accurate pointwise approximation of the unknown function to be 
recovered on certain segments of the domain Q. We next illustrate on a variety of PDE 
problems and OCPs, that indeed good pointwise approximation can be obtained on some 
parts of the domain Q. 



5 Numerical Experiments 



We illustrate the potential of Method [2] on a variety of ODE and PDE problems. As an imple- 
mentation of the sparse SDP relaxations we apply the software SparsePOP of Waki et al.(2005)] 
and as an implementation of sequential quadratic programming (SQP) in the SDPR method 
and in solving the optimization problem ([7]) we apply the Matlab Optimization Toolbox com- 
mands fmincon and fminunc, respectively. As we restrict ourselves to ODEs and PDEs 
with rectangular domains fl we can apply standard quadrature and cubature formulas to 
compute the gradient and Hessian for ([7]). Thus, we apply the Matlab commands trapz 
in the one- dimensional case and dblquad in the two-dimensional case, respectively. In or- 
der to evaluate the quality of the smooth approximation provided by Method [2] we de- 
fine the average error eu{M) := AxJ^i'^i^i) ~ '^X-ii.^i) ^-^^ the maximum pointwise error 
gmax^^^ := maxi | u{xi) - ulj{xi) |. 



5.1 Linear differential equations 

As first test problems for our approach we consider a linear ODE and a linear PDE from 
Bertsimas et al. (2006)| and compare (a) the discretized moment approximations obtained 
by Method [2] to the bounds obtained in [Bertsimas et al. (2006)| , and (b) the smooth approx- 
imation u*^ to the known analytic solution. 



5.1.1 Linear ODE 



The linear ODE is given by, 

u"{x) + 3u'{x) + 2u{x) = Va;e[0,l] 



u'{0) = -2e\ u'{l) = -2. 

For this problem the unique solution is u{x) = e^e'"^^ . We apply Method [T] with F = — J2i ^i, 
N = 2000 and w = 1. We calculate the approximate moments for M = 40 and compare 
them to approximate moments for m derived in Bertsimas et al. (2006)] by contracting lower 



and upper bounds in Table [T] All moments coincide up to the fourth digit. 

We apply Method [2] for M G {1, . . . , 5} and report the resulting vector v^,j in Table |2J 

Note that the actual solution of the ODE problem corresponds to a parameter vector fopt = 
(2, -2, 0, . . .). As reported in Table El {v^{M),v^{M)) ^ (2, -2), and even though v*{M) ^ 
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'"'t 


'"'i 


n 


3.1942 


3 1 Q45 


1 




1 0Q73 


10 


0.1086 


0.1088 


20 


0.0524 


0.0524 


30 


0.0345 


0.0344 


40 


0.0258 


0.0256 



Table 1: Approx. moments for = 2000 compared to approx. m^*" derived from lower 



and upper bounds in Bertsimas et al. (2006)| for linear ODE 



M 


1 


2 


3 


4 


5 


vl 


2.0043 


2.0042 


2.0040 


2.0032 


2.0039 


v\ 


-2.0135 


-2.0083 


-1.9922 


-1.9744 


-1.9528 






-0.0083 


-0.0938 


-0.1998 


-0.3205 


vl 






0.0851 


0.2886 


0.3269 


vl 








-0.1180 


0.2018 


vl 










-0.2694 


euiM) 


0.0042 


0.0044 


0.0042 


0.0039 


0.0044 




0.0242 


0.0228 


0.0228 


0.0222 


0.0326 



Table 2: Smooth SDP approx. for linear ODE (Q 



for i > 1, the maximum pointwise error is quite small. However, {vl{M),vl{M)) does not 
converge to (2, —2) since pointwise convergence is not guaranteed. In fact the higher moment 
terms v*{M) for i > 1 counterbalancce the difference between {vl{M),vl{M)) and (2, —2). 



5.1.2 Linear PDE 



The linear PDE is given by 

u.x{x, y) + Uyy{x, y) - 3e"+^ = V (a;, G [0, 1]^, (10) 

where the boundary conditions are set up such that u{x^ y) = e^"*"^ is the unique solution 
of the PDE problem. We apply the SDPR method with F{u) = —Ylk I'^k,!, 1^1 = 1 and 
Nx = Ny = 100, compute for all moments of order up to 3 and compare these moments to 
the exact moments derived from the known solution and to the best lower and upper bounds 
for the moments from Table 7 in Bertsimas et al. (2006)| . See Table [3] for the results. Our 



approach provides some approximations for the moments that is well within the bounds LB 
and UB from Bertsimas et al. (2006)| . 



Also for this problem, there is a vector t>opt = (0, 1, 1, 0, . . .) corresponding to the solution u. 
When applying Method [2] with M G {2,3} we obtain v*{M) reported in Table H] with errors 
e„(2) = 0.0677 , e„(3) = 0.0682, ej^^^(2) = 1.0276 and ej^^^(3) = 1.3267. 

For both linear differential equations we obtain accurate approximations of the moments for 
the measures associated with the unique solutions. Moreover, unlike Bertsimas et al. (2006)| , 
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V*") J ) 




771 • ■ 


LB 


UB 


fn {]) 


2.9512 


2 9525 


2 9235 


3 1 707 




1 7287 


1 7183 


1 6944 


1 7742 

X.I 1 




1.0121 


1.0000 


0.9847 


1.0130 


(2,0) 


1.2504 


1.2342 


1.1123 


1.3088 




7319 


71 83 


7151 

W.I X X 


7458 


(2 2) 


5292 


51 59 


4948 


5456 


(3,0) 


9869 


0.9681 


0.8244 


1.0478 


(3,1) 


0.5775 


0.5634 


0.5054 


0.5818 


(3,2) 


0.4176 


0.4047 


0.3874 


0.5818 


(3,3) 


0.3295 


0.3175 


0.3017 


0.3399 



Table 3: Approx. moments for A^^^^ = Ny = 100 compared to exact moments and bounds from 
Bertsimas et al. (2006)] for linear PDE ^ 





(0,0) (1,0) (0,1) (2,0) (1,1) (0,2) 


<,(2) 


0.064 0.709 0.710 0.329 0.099 0.329 


<,(3) 


0.013 0.959 0.959 -0.213 0.160 -0.213 




Oil 00 



Table 4: Smooth SDP approx. for linear PDE ( fTOj) 



we can exploit these moment approximation to find smooth approximations for the actual 
solutions. 



5.2 A nonlinear elliptic PDE 

As a first nonlinear problem, consider the elliptic PDE: 

Uxx + Uyy + 22w(l -u^) = on [0, l]^ 

u = on9[0, 1]2, (11) 

< M < 1 on [0, 1]2. 

This problem is well known to have a nontrivial positive solution. We apply Method [2] with 
F{u) := — ^kiUk,u w = 2, Nx = Ny = 49 and M G {2,3}, and obtain an accurate discrete 
approximation {uij)ij and two smooth approximations U2 and with errors (e™^^(2), 6^(2)) = 
(0.081, -0.002) and (e5^''=^(3), e„(3)) = (0.076, -0.002), respectively. In FigureODit is illustrated 
how accurate the smooth approximation resembles the shape of the discrete approximation 
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Figure 1: Discrete approx. (left), smooth approx. (center) for the same grid, and pointwise 
difference between the two (right, notice the vertical axis scale) for nonlinear PDE (1 11 1) 

5.3 Reaction diffusion equation 



A challenging ODE problem which is known to have many solutions is given as in |Mimura(1979)] 
by: 



u" + 1(35 + 16u - 1/2) u-uv = ^o\i [0,5], 
iv" - {\^\v) V ^uv = ^ on[0,5], 
v'(0) = v'(5) = 0, 



(12) 



0<M,v<14 on [0,5]. 

To ensure numerical stability we scale the domain [0, 5] of (IT^ to [0, 1] before applying Method 
El as exp(x*) gets very large for > 1. Problem f lT2|) involves two functions u and f, thus 
we need to apply Method [2] twice, once for u and once for f , in order to obtain smooth 
approximations for both functions. Note that the SDPR method only needs to be applied 
once and provides a discrete approximation {ui,Vi)i, maximum entropy estimation needs to 
be applied twice. We apply Method |2] for F{u,v) = —u^n^ , = 100, w = 3 and M G 
{10, 20, 30, 50}. Solutions for u and v and the derived smooth approximations are pictured in 
Figure [2] and the errors are reported in Table [51 Even though the average and the maximum 
pointwise approximation errors are not small, even for large M, we observe that the pointwise 
approximation gets more and more accurate on certain segments of the domain - in this 
case the third part of the interior, where the third peak of u and v occurs. The pointwise 
approximation is not good - or does even get worse - near the boundary of [0, 1]. 





Figure 2: Smooth approx. for u (left) and v (right) for M G {10,20,30,50} for nonlinear 
ODE ([H 
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IVl 


_max/ i\/r\ 






ey{ivi) 


1 n 


4 1 1 


n DD1 1 


9 1 n 


n 1 DRD 

U. lUOU 


20 


3.84 


0.0143 


6.66 


0.1371 


30 


3.48 


0.0012 


9.61 


0.1547 


50 


3.39 


0.0085 


23.69 


0.2658 



Table 5: Errors for nonlinear ODE flT^ 



5.4 Control of production and consumption 



A first example of an optimal control problem is given by 

rT 



min — Jq (1 — u{t))x{t) dt 
s.t. x{t) = u{t)x{t) Vt e [0, T] , 

x(0) = xo, (13) 
< x{t) < 1 Vt G [0, 1] , 

< xlt) < 10 Vt G [0, 10] , 



where T > 1 fixed. For this simple problem the optimal control law is given by 



1 ifO<t<T-l, 
if r - 1 < t < T. 



We choose xq = 0.25, T = 4 and = 100, and apply Method [2] after scaling the domain to 
[0, 1] for M G {5, 10, 20, 30}. In this case the objective function F is given by a discretization 
of the objective function of f|T3l) . As in the case of ordinary and partial differential equations, 
we observe that a fairly good pointwise approximations of both, optimal control law and 
corresponding trajectory, are obtained on the interior of the domain. See Figure [31 




11 



5.5 The double integrator 



Another interesting control problem, which has been discussed in Lasserre et al. (2008)1 , IS 
given by 



mm 
s.t. 



T 

X2it) 

x{0) -- 
x{T) 



-- X2{t) 

--u{t) 
= (0,0), 



-1 < u{t) < 1 

-1 < Xi{t),X2it) < 10 



yt E [0,T] 

Wt e [o,T] 



vtG [o,r], 
vt G [o,r]. 



(14) 



After scaling the domain to [0, 1], we apply Method|5]with = 50, a; = 3, Xq = (0.8, —1) and 
M e {5, 10, 20, 30}. Since T is not fixed, it is treated as an additional variable, i.e., we consider 
the polynomial optimization problem with variable {xn, ... ,xin,X2i, ■■■ ,X2n,ui, ... ,un,T). 
We observe that the pointwise approximation of optimal control and trajectories on the interior 
of the domain gets better and better for increasing moment order, see Figure |H Note, unlike 
the moment based method in Lasserre et al.(2008)| which yields bounds for the optimal value 
of the control problem, we obtain discrete and smooth approximations to the optimal control 
and the corresponding trajectories by Method [1] and Method El respectively. 




Figure 4: Smooth approx. for Xi (left), X2 (center) and u (right) for M G {5,10,20,30} 
forOCP (dl. 



6 Concluding Remarks 



We introduced a novel technique to derive smooth approximations for solutions of systems of 
differential equations and optimal control problems, which is based on sparse SDP relaxations 
and the maximum entropy estimation method. As demonstrated on some examples of nonlin- 
ear ordinary differential equations, partial differential equations and optimal control problems, 
this technique is promising to obtain accurate pointwise approximations of a solution of a dif- 
ferential equation on segments of its domain. It would be interesting to characterize regions of 
the domain of nonlinear differential equations where accurate pointwise approximations can 
be guaranteed. Another question is whether a different, less restrictive parametrization than 
(|6]) allows better pointwise approximations on larger segments of the domain. Of course, other 
choices of entropy would lead to different results, and it remains a topic of future investigation 
to analyze, whether or not some choice should be preferred. At this point solving the sparse 
SDP relaxation is the major computational bottleneck in Method [2j Every improvement for 
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approximation accuracy and efficiency of SDP relaxation techniques for polynomial programs 
and every improvement for SDP solvers will extend the range of the proposed method. 

In order to validate these techniques in the present context of PDEs, a comparison with state- 
of-the-art numerical methods for solving nonlinear PDEs remains to be done. In particular, 
the maximum entropy technique should be compared with standard and spline interpolation 
methods which both provide differentiable approximations that coincide with the discrete so- 
lution at each grid point. But, unlike our method, in polynomial interpolation the degree 
of the polynomial is directly related to the number of grid points. Also, cubic splines are 
twice differentiable at the grid points, whereas our approximation is smooth. In contrast, and 
instead of searching for a smooth solution that matches the unknown u at the grid points, 
one searches for a smooth solution that matches finitely many moments of the associated 
measure = udx. And so, an interesting feature is to obtain a smooth approximation from 
a limited number of moments and not from the discrete solution explicitly (still, before using 
moments, a sufficiently fine discretization is necessary to ensure ^ m). This method is an 
attempt to deal with the curse of dimensionality in the numerical analysis of PDEs. Further- 
more, we emphasize that the proposed method can be applied to a wide range of problems 
that encompasses nonlinear partial differential equations and nonlinear optimal control prob- 
lems. A detailed comparison is however out of the scope of this paper, whose main objective 
was to pave the way for the development of SDP techniques to solve optimal control design 
problems for dynamical systems described by polynomial ODEs and PDEs, in the spirit of 
Lasserre et al.(2008)| . 
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